#ifndef	_CONTAINMENT_H_
#define _CONTAINMENT_H_

#include <numeric>

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Aff_transformation_2.h>


namespace Packing {

	class Constraint
	{
		typedef CGAL::Exact_predicates_inexact_constructions_kernel K;

	public:
		typedef K::Point_2 Point_2;
		typedef K::Segment_2 Segment_2;
		typedef K::Vector_2 Vector_2;
		typedef CGAL::Aff_transformation_2<K> Transformation_2;

	public:
		Constraint(const Point_2& p, const Segment_2& seg, const Point_2& ref_point): p_(p), dist_(0.0), seg_(seg)
		{
			Vector_2 ref_v(seg_.source(), ref_point);
			n_ = seg_.to_vector().perpendicular(CGAL::COUNTERCLOCKWISE);
			if (n_*ref_v < 0)
				n_ = -n_;
			n_ = n_ / CGAL::sqrt(n_.squared_length());
		}
		Constraint(const Point_2& p,const Segment_2& seg ) : p_(p), dist_(0.0), seg_(seg)
		{
			Vector_2 v(seg_.source(), p_);
			n_ = seg_.to_vector().perpendicular(CGAL::COUNTERCLOCKWISE);
			if ( v*n_ < 0)
				n_ = -n_;
			n_ = n_/sqrt(n_.squared_length());  // normalize
		}
		Constraint(const Point_2& p,const Segment_2& seg, double dist ): p_(p), dist_(dist), seg_(seg)
		{
			//assume p has been already inside the region when initialization
			Vector_2 v(seg_.source(), p_);
			n_ = seg_.to_vector().perpendicular(CGAL::COUNTERCLOCKWISE);
			if ( v*n_ < 0)
				n_ = -n_;		
			n_ = n_/sqrt(n_.squared_length());  // normalize
		}
	public:
		Point_2 p_; // the transformed p should be on the positive side of the segment
		Segment_2 seg_;
		double dist_;
		Vector_2 n_; // unit vector perpendicular to the segment

		double transformed_signed_area(const Transformation_2& aff)
		{
			Point_2 pp = aff(p_);
			return (pp - seg_.source())*n_;
		}
		double transformed_signed_dist(const Transformation_2& aff, double k)
		{
			Point_2 pp = aff(p_);
			return (pp - seg_.source())*n_ - k*dist_; 
		}
		void transform(const Transformation_2& aff)
		{
			p_ = aff(p_);
			seg_ = seg_.transform(aff);
		}
	};

	class Containment
	{
		typedef Constraint::Point_2 Point_2;
		typedef Constraint::Segment_2 Segment_2;
		typedef Constraint::Vector_2 Vector_2;
		typedef Constraint::Transformation_2 Transformation_2;

	public:

		void insert_constraint(const Point_2& p,const Segment_2& seg, double dist);  // the distance from the transformed p to the line should be greater than dist

		void compute_constraints(const double * const x, double* c, int mc);  // interface for Knitro

		void compute_constraint_grads(const double * const x, double* jac);

		void centralize_constraints(const Point_2& local_origin);

		int get_constraint_list_size(){return const_list_.size();}

		void clear();

	public:
		std::vector<Constraint> const_list_;

	};


}

#endif